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ABSTRACT 


This  report  discusses  the  implementation  of  one  of  the  best  digital  cyclic  spec¬ 
trum  (CS)  algorithms  derived  so  far,  the  Strip  Spectral  Correlation  Algorithm  (SSCA). 
Some  theoretical  background  and  a  detailed  description  of  the  SSCA  are  provided.  An 
analysis  of  the  SSCA  is  performed  and  an  algorithm  for  mapping  the  SSCA  output  is 
formulated.  The  cyclic  feature  function  (CFF)  is  defined  as  a  means  to  detect  the  cyclic 
features  from  the  SSCA.  Results  of  the  SSCA  encoded  in  C  are  then  reported.  Three 
BPSK  signals  and  two  additive  white  Gaussian  noise  (AWGN)  signals  are  used  to  verify 
the  validity  of  the  SSCA.  Three-dimensional  plots  and  two-dimensional  plots  of  the  CS 
and  CFF  respectively  are  presented  to  the  reader.  Finally,  some  benchmarks  on  a  SUN 
computer  for  the  SSCA  are  provided  for  reference.  In  brief,  the  CS  and  CFF  estimated 
with  the  SSCA  prove  to  be  valuable  tools  for  analyzing  second-order  cyclostationary  com¬ 
munication  signals  and,  by  making  extensive  use  of  the  FFT,  to  provide  robust,  reliable, 
and  accurate  results  more  efficiently  than  typical  CS  direct  estimation  methods. 

RESUME 

Ce  rapport  traite  de  la  mise  en  oeuvre  de  I’un  des  meilleurs  algorithmes  elabores  a 
ce  jours  pour  ce  qui  est  de  I’estimation  complete  du  spectre  cyclique  (SC)  d’une  sequence 
numerique,  soit  I’algorithme  SSCA.  Les  bases  theoriques  et  la  description  detaillee  du 
SSCA  sont  fournies.  L’analyse  du  SSCA  est  effectuee  et  un  algorithme  destine  a  reordonner 
les  resultats  du  SSCA  est  formule.  La  function  de  caracteristique  cyclique  (CFF)  est 
definie  et  employee  afin  de  detecter  les  caracteristiques  cy cliques  a  partir  de  I’algorithme 
SSCA.  Les  resultats  du  SSCA  code  en  langage  C  sont  ensuite  rapportes.  Trois  signaux 
BPSK  ainsi  que  deux  signaux  composes  de  bruit  blanc  Gaussien  additif  sont  utilises  pour 
verifier  la  validite  du  SSCA.  Les  representations  graphiques  en  deux  et  trois  dimensions 
respectivement  du  CFF  et  du  SC  sont  presentees  au  lecteur.  Finalement,  des  points  de 
reperes  en  ce  qui  concerne  le  temps  d’execution  du  SSCA  sur  un  ordinateur  SUN  pour 
les  parametres  utilises  sont  fournis  a  titre  de  reference.  En  bref,  le  SC  et  le  CFF  estimes 
a  I’aide  du  SSCA  constitue  assurement  des  outils  importants  pour  I’analyse  de  signaux 
cyclostationnaires  de  deuxieme  ordre  et,  par  une  vaste  utilisation  de  la  FFT,  le  SSCA  per- 
met  de  produire  des  resultats  robustes,  fiables  et  precis  de  fagon  beaucoup  plus  efficace  ; 
que  d’autres  methodes  numeriques  plus  directes  d’estimation  du  SC. 


EXECUTIVE  SUMMARY 


Athough  much  superior  to  ordinary  spectral  analysis,  cyclic  spectral  analysis 
requires  considerably  more  complex  and  time-consuming  computations  to  produce  an  es¬ 
timate  of  the  entire  cyclic  spectrum  (CS).  An  estimate  of  the  entire  CS  is  required  when 
a  priori  information  on  signal  cycle  frequencies  is  not  available.  In  this  case,  detection 
of  cycle  frequencies  possibly  through  the  CS  is  crucial  to  any  cyclostationary  property 
exploiting  methods  for  several  applications  encountered  in  typical  military  surveillance 
system  such  as,  for  example,  signal  detection  and  classification,  signal’s  parameter  esti¬ 
mation,  time-difference  of  arrival  estimation,  adaptive  spatial  filtering,  direction  finding 
and  signal  extraction.  Based  on  either  a  frequency  or  time-smoothing  approach,  efficient 
algorithms  for  estimating  the  cyclic  spectrum  of  discrete-time  domain  signals  have  been 
devised  during  the  last  decade  or  so.  This  report  discusses  the  implementation  of  one  of 
the  best  digital  CS  estimation  algorithms  derived  so  far,  the  Strip  Spectral  Correlation 
Algorithm  (SSCA). 

Following  a  review  on  the  theory  of  spectral  correlation  for  cyclic  spectral  analy¬ 
sis  of  discrete-time  signals,  the  SSCA  is  derived.  The  detailed  description  shows  that  the 
SSCA  essentially  reduces  in  complexity  to  a  two-dimensional  FFT  of  size  NxNp  and  that 
it  is  suitable  for  parallel  architectures.  Parameters  N,  denoting  the  observation  time,  and 
Ap,  specifying  the  number  of  samples  considered  in  the  channelizer,  determine  a  number  of 
CS  estimate  attributes  such  as  the  frequency  resolution  A/,  the  cycle  frequency  resolution 
Aa,  and  the  time-frequency  resolution  product  A/At.  All  the  steps  required  to  imple¬ 
ment  the  SSCA  are  clearly  explained.  An  analysis  of  the  SSCA  reveals  its  computional 
requirements.  The  accuracy  of  the  estimates  are  also  considered  and  an  algorithm  for 
mapping  the  SSCA  output  is  provided.  Additionally,  the  cyclic  feature  function  (CFF) 
is  defined  as  a  means  to  detect  the  cyclic  features  from  the  SSCA.  Results  of  the  SSCA 
encoded  in  C  are  then  reported. 

Theoretical  CSs  axe  derived  for  three  BPSK  signals  using  various  parameters  and 
two  additive  white  Gaussian  noise  (AWGN)  signals.  Three-dimensional  plots  of  the  CS 
and  two-dimensional  plots  of  the  CFF  of  various  combinations  of  these  signals  along  with 
various  algorithmic  parameters  are  realized  with  the  ssca  program  and  presented  to  the 
reader.  Specifically,  the  effect  of  varying  the  observation  time  and  the  frequency  resolution 
are  observed  on  one  BPSK  signal  and  one  AWGN  independently.  Then,  the  CS  and  CFF 
of  combinations  of  the  signals  under  investigations  are  studied  while  keeping  all  SSCA 
parameters  fixed.  Cyclic  features  are  detected  based  on  the  CFF.  These  experimental 


V 


results  show  that  all  major  cyclic  features,  corresponding  directly  or  indirectly  to  signal 
features  (e.g.  carrier  frequency  and  baud  rate  in  this  case),  are  successfully  detected 
even  in  highly  corrupted  conditions  where  conventional  analysis  would  have  failed  to 
discriminate  the  signals.  Finally,  some  bechmarks  on  the  SUN  computer  for  the  SSCA 
are  provided  for  reference. 

In  brief,  the  CS  and  CFF  estimated  with  the  SSCA  prove  to  be  valuable  tools  for 
analyzing  second-order  cyclostationary  communication  signals  and,  by  making  extensive 
use  of  the  FFT,  to  provide  robust,  reliable,  and  accurate  results  more  efficiently  than 
typical  CS  direct  estimation  methods.  This  suggests  further  work  to  apply  the  SSCA 
to  off-air  data  and  possibly  implement  the  algorithm  in  a  popular  digital  signal  multi¬ 
processor  board  such  as  a  C40  board. 
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1.0  INTRODUCTION 


Cyclic  spectral  analysis  (CSA)  as  a  general  tool  for  signal  analysis  has  been  shown 
to  be  much  superior  to  more  conventional  spectral  analysis  [1].  However,  it  is  recognized 
that  the  computation  of  the  entire  cyclic  spectrum  (CS)  is  considerably  more  complex  and 
time-consuming  than  the  computation  of  the  conventional  power  spectral  density  (PSD). 
The  estimation  of  correlation  between  spectral  components  of  signals  as  compared  to 
only  computing  the  spectral  components  themselves  makes  the  CSA  a  computationally 
complex  mechanism  mainly  because  of  the  potentially  large  number  of  correlation  com¬ 
putations  involved.  With  the  growing  importance  of  CSA  for  many  applications  (cf.[2]), 
a  search  to  derive  efficient  algorithms  for  computing  the  cyclic  spectrum  has  been  taking 
place  over  the  last  decade  or  so.  These  CS  estimation  algorithms  divide  mainly  into  two 
subgroups:  the  frequency-smoothing  (FS)  algorithms  and  the  time-smoothing  (TS)  algo¬ 
rithms.  More  lately,  algorithms  merging  the  two  methods  have  been  used  and  are  referred 
to  as  hybrid-smoothing  algorithms  [3].  These  will  not  be  discussed  in  this  report. 

Focusing  on  the  problem  of  computing  the  entire  CS,  TS  algorithms  result  in 
more  computationally  efficient  algorithms  than  those  involving  FS.^  Based  on  an  analysis 
of  the  performance  of  the  currently  existing  TS  algorithms  (in  [4]  and  [3]),  one  particular 
TS  algorithm  called  the  Strip  Spectral  Correlation  Algorithm  (SSCA)  has  been  chosen 
for  implementation. 

This  report  looks  at  the  implementation  of  the  SSCA  for  the  computation  of 
the  entire  CS.  But  first,  some  background  on  CS  estimation  for  discrete-time  signals  is 
required. 


2.0  BACKGROUND 


In  this  section,  the  theory  of  spectral  correlation  necessary  for  CSA  of  discrete¬ 
time  signals  is  reviewed.  Fundamentally,  there  are  two  possible  approaches  for  developing 
the  theory  of  CSA:  the  probabilistic  approach  (cf.[5])  and  the  deterministic  approach 
(cf.[6]).^  Here,  the  deterministic  framework  is  chosen  because  of  its  greater  tractability 

^Note  that  if  the  CS  is  required  for  only  few  values  of  known  cycle  frequency,  FS  algorithms  can  be 
computationally  superior. 

^In  the  deterministic  framework,  time  averages  are  used  as  opposed  to  ensemble  averages  for  the 
probabilistic  framework. 


and  its  natural  link  to  practical  problems  [2]  (e.g.,  signal  processing  of  discrete-time 
communication  signals  for  various  applications). 


The  cyclic  autocorrelation  function  of  a  discrete-time  sequence  a;(n)  is  defined 


The  sequence  x{n)  is  said  to  be  wide  sense  cyclostationary  if  and  only  if  R^ik)  ^  0  for 
some  Q.  The  CS  S^{f)  is  simply  the  Fourier  Series  Transform  (FST)  of  R°{k),  i.e., 

StU)  =  E  .  (2) 

k=—oo 

Defining  a  as  the  cycle  frequency  and  /  as  the  spectral  frequency,  S^{f)  can 
give  rise  to  a  three-dimensional  plot  with  the  three  axes  being:  either  |<S’“(/)|  or  Z5"(/), 
/,  and  a.  In  general,  S^(f)  is  complex- valued. 

2.1  Basic  Frequency-smoothing  Approach 

The  FS  approach  for  computing  S°{f)  is  based  on  the  double  limit 

=  (3) 


where 


S;jnJU,  =  XMn,F  +  a/2)Xi,(n,F-a/2)dF 


is  the  frequency-smoothed  cyclic  periodogram  and 


m^n  2 


is  referred  to  as  the  complex  demodulate  of  x(n),  a  discrete-time  signal  sampled  at 
/.  =  1/Z  Hz. 

5'"^^(n,/)A/  is  the  result  of  correlating  spectral  components  of  x{n)  over  A/  Hz 
(i.e.,  averaging  in  frequency).  The  quantities  A/  and  At  are  respectively  called  the  fre¬ 
quency  and  time  resolution  of  the  estimate.  Note  that  f)^/  is  used  for  estimating 
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S°{f)  and  a  reliable  estimate  would  have  a  large  resolution  product  AtAf. 


2.2  Basic  Time-smoothing  Approach 


Using  the  TS  approach,  5"(/)  is  also  expressed  as  a  double  limit 


i  — ^OO  At— ^00  ■* 


where 


=  E  XTim,f  +  al2)X*AmJ-al2) 


(6) 


(7) 


is  the  time-smoothed  cyclic  periodogram  and 

1  n+f 

XT{nJ)  =  7f  E 


(8) 


is  another  form  of  the  complex  demodulate  of  x{n)  for  a  smaller  time  duration  since 
T  <  At. 


The  TS  cyclic  periodogram  is  used  to  estimate  <S'"(/).  As  can  be 

seen  from  (7),  S“^{n,  f)^t  is,  in  fact,  the  discrete  time-average  of  spectral  correlation 
components  of  x(n)  over  a  time  At.  In  this  case,  the  frequency  resolution  is  Af  =  1/T. 


2.3  Comments  on  Both  Approaches 

Some  comments  with  respect  to  both  approaches  have  to  be  made. 

1.  The  order  of  the  double  limit  used  to  compute  5'“(/)  is  not  interchangeable. 

2.  It  has  been  shown  in  [6]  that  FS  and  TS  approaches  are  equivalent  if  and  only  if 
AtAf  >  1,  i.e.,  5“^,(n,/)A/  «  S°^{n,  f)^t  for  AtAf  >  1. 

3.  The  condition  AtAf  >•  1  is  essential  in  order  to  obtain  reliable  estimates. 

4.  The  particular  choices  of  AtAf  and  Af  depend  on  the  signals  under  investigation 
(signals  environment)  and  the  degree  of  reliability  which  is  needed. 


3 


Generally  speaking,  the  goal  is  to  find  an  efficient  algorithm  for  high  At  A/.  In 
[4],  it  is  shown  that  the  SSCA  (a  TS  algorithm)  is  one  of  the  best  derived  so  far  for 
AtAf  >  128,  Af  >  1/32,  and  with  a  small  Aa  (the  cycle  frequency  resolution). 

3.0  THE  STRIP  SPECTRAL  CORRELATION  ALGORITHM 

This  section  describes  and  analyses  the  SSCA. 

3.1  From  the  Basic  Time-smoothing  Approach  to  the  SSCA 

Using  a  data  tapering  window  a(r)  of  length  T  —  NpTs  seconds,  the  complex 
demodulate  of  x(n)  can  be  re-expressed  as  being 

XT(nJ)  =  E  o(r)i(n  - 
—  ^  a{r)x{n  + 

-Np 

Lr=-5^ 

where  T,  is  the  sampling  period.  The  term  in  brackets  represents  the  centered  Fast  Fourier 
Transform  (FFT)  of  x(n)  using  a  window  a(r).  Note  that  Af  =  Aa  =  fsiNp  (where 
Aa  =  5:  is  the  bandwidth  of  the  input  filters).  Correlating  the  complex  demodulates  for 
a  time  At  =  NT3  seconds,  the  TS  cyclic  periodogram  is  then  given  by 

S°^{nJ)At=  X)  (n  +  m,  f  +  Xt  (n  +  mj  -  g{m)  (11) 

m=-f 

where  g{m)  is  another  data  tapering  window  of  length  At.  As  mentioned  earlier,  for 
S^^{n,f)At  to  be  a  reliable  estimate  of  AtAf  >■  1  (i.e..  At  T)  is  required.  In 

the  limit  (i.e..  At  — >  00  followed  by  Af  0),  if  a(n)  and  g{n)  are  normalized  such  that 

X  =  Zi  9{n)  =  1 ,  (12) 

n  n 

then  5"^(n,/)At  =  5“(/)  [3]. 

In  order  to  extend  (11)  to  transform  the  averaging  operation  into  an  FFT  ,  the 


r)T,  (9) 
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product  sequence  must  be  frequency  shifted  by  an  amount  7  (i.e.,  from  a  to  a  +  7).  The 
result  is  then 


^-1 

SJ+>,/U,=  E  +  +  (13) 

Using  (13),  one  can  find  that  Af  =  Aa  —  \'y\  with  (I7I  <  Ao),  a  non-uniform  frequency 
resolution  since  several  values  of  7  are  usually  required.  If  7  =  qAa,  then  (13)  becomes 


—j27rmq 


5’“;+’^“(n,  fj)At  =  2  XT(n  +  m,  fk)  Xt  (n  +  m,  //)  ^(m)e  ^ 


Tn=— 


which  is  an  N-point  centered  FFT.  Note  that  fk  =  k{fsfNp),  k  =  —Np/2,  (A'p/2)  —  1, 

represents  several  frequencies  for  a  complete  coverage  of  the  bifrequency  plane  (/,  a).  To 
relate  fj  and  a,-  to  fk  and  //,  the  following  relations  are  needed: 

r  fk  +  fl  k  +  I  f  fs\  ,  , 


2  \K 


For  the  Strip  Spectral  Correlation  Algorithm  (SSCA)  [3][4][7][8],  one  has  to  set 
Q,-  =  fi^^  fj  =  (^fk  —  qAa)  12  and  replace  Xj{n,fi)  by  the  unfiltered  signal  x*(n)  in  (14)^ 
such  that 

— -1 

^fk+9^‘^  L  ^  XT(n  +  m,fk)x*(n  +  m)g(m)e-^^.  (17) 

\  ^  /At  rr,=-f 

This  allows  Af  as  well  as  AtAf  to  be  uniform,  a  highly  desirable  feature  for  the  estima¬ 
tor.  The  points  estimates  produced  by  (17)  will  lie  along  the  frequency-skewed  family  of 
lines  a  =  2fk  —  2f  to  result  in  a  strip  for  each  /*.  The  following  applies  to  the  SSCA 


^Note  that  this  is  done  at  the  expense  of  minor  degradation  in  output  signal-to-noise  ratio  (cf.[7]). 
The  estimate  obtained  in  (17)  is  equivalent  to  (14)  only  if  A/At  >  1,  condition  which  is  essential  to  have 
a  reliable  estimate  anyway. 


implementation: 


A/  =  l/7Vp,  Aa  =  1/At  =  1/N,  AfAt  =  N/Np,  a(0)  =  XI S'!”)  =  (1^) 

n 

and 

/fc  +  ^Aa,  k  =  -^,  ...,^-l  (19) 

fk  —  qAa  N  N 

Note  that  if  the  cross-cyclic  spectrum  (say  between  x(n)  and  y{n))  is  required  instead 
of  the  auto-CS,  then  x*{n  m)  in  (17)  simply  becomes  y*{n  +  m).  In  addition,  if  the 
conjuguate  auto-CS  and  conjuguate  cross-CS  are  desired,  simply  replace  x*{n  -|-  m)  and 
y*{n  -f  m)  by  x{n  -}-  m)  and  y{n  -f  m)  in  (17)  respectively.  It  has  been  pointed  out  in  [9] 
that  since  one  of  the  channels  is  unfiltered,  the  complex  demodulates  have  to  come  in  at 
full  rate,  i.e.  without  decimation.  It  is,  however,  proposed  in  that  particular  paper  to 
decimate  the  complex  demodulates  in  order  to  reduce  the  number  of  input  channelizers 
by  incorporating  a  hold  operation  at  the  end  of  every  channelizer.  This  method  has  been 
tried  here  by  the  author  and  it  was  found  that  the  hold  operation  was  introducing  false 
cyclic  features  which,  on  the  other  hand,  could  possibly  be  predicted.  The  technique 
will  not  be  discussed  further  in  this  paper  but  merits  some  considerations  as  it  saves 
computation  time  and  storage  space. 


3.2  SSCA  Description 

Based  on  (10)  and  (17)  and  assuming  a  single  processor,  the  algorithm  can  be 
formulated  in  the  following  steps. 

Step  1  Data  collection: 

•  Collect  a  block  of  data  of  iV  -|-  iVp  samples,  i.e.,  x(n),  n  =  1,  2,  . . . ,  -f  iVp  . 
Step  2  Compute  N  jV„-point  FFTs  of  sub-blocks  of  x(to): 

•  Define  a{r)  for  r  =  1,  2,  . . . ,  iVp  to  be  an  iVp-point  data-tapering  window  (e.g. 
rectangular,  Hamming,  Hanning,  Kaiser  ...). 


6 


•  For  n  =  1,  2,  . . . ,  iV,  r  =  1,  2,  . . . ,  Np,  and  k  =  — . . . ,  ^  —  1,  compute 

XT{n,  k)  =  FFTStVp  {a{r)x{n  +  r  —  1)}  (21) 

where  FFTS^Vp  is  the  iVp-point  centered  FFT  operation  (i.e.,  with  the  zero 
frequency  in  the  k  =  0  bin). 

Step  3  Compute  the  weighted  product  from  (17): 

•  Define  g(n)  for  n  =  1,  2,  . . . ,  iV  to  be  a  A^-point  data-tapering  window  (e.g. 
as  above). 

•  For  n  =  1,  2,  . . . ,  N  and  k  =  — •  - . ,  ^  —  1,  compute 

corrected  FFTs  ^  in  (10) 


Xg{n,  k)  =  XT{n,  k)e  ^  x*  j 


Step  4  Compute  the  spectral  correlation  function: 

•  For  q  =  — y,  . . . ,  y  —  1,  fc  =  — . . . ,  ^  —  1,  and  n  =  1,  2,  . . . ,  compute 

S,{q,k)  =  FFTSN{Xg{n,k)}  (23) 

where  FFTSjv  is  the  A^-point  centered  FFT  operation  (i.e.,  with  the  zero  fre¬ 
quency  in  the  q  =  0  bin). 

Step  5  Map  5j(<y,  k)  onto  S°(f): 

•  The  mapping^  is  performed  using  the  equations 


f  —  -h _ 2_ 

J  ~  2Np  2N 

a  =  —  4-  -2- 

^  —  AT.  I  ?<I 


where  /  and  a  are  normalized  with  respect  to  /,  =  1,  i.e., 


0.5  </<  0.5  and  -  1  <  a  <  1 


^This  mapping  is  absolutely  required  for  plotting  and  interpreting  the  cyclic  spectrum. 

®Note  that,  really,  fmin  ^  f  ^  fmax  where  fmin  “  0.5  -b  2/7  and  fmaz  —  O.b  2/v  while 


O^min  ^  ^  ^  ^max  where  OCfnin  —  1  and  (Xfnax  —  1 


27V,  2N  • 


7 


This  algorithm  (as  described  from  Step  1—5)  can  also  be  expressed  in  terms  of 
matrix  operations  (or  array  processing).  One  possible  formulation  is  presented  below. 

Step  1  Given  the  input  data  vector  x=  a;(l)  x{2)  •••  x{N  +  Np)  ],  create  the 
NxNp  input  data  matrix  X  where 


a;(l)  x{2) 

x(2)  a;(3) 


x{Np) 
x{Np  +  l) 


[x(iV)  x{N+l)  •••  x{N  +  Np)\ 

Step  2  Compute  the  TVp-point  window  row-vector  a=  |  a(l)  a(2)  •••  a{Np)  j  and 
create  A  (an  NxNp  matrix)  formed  by  N  rows  of  a.  Knowing  A,  compute  XA  =  X  ©  A 
where  0  denotes  the  Hadamard  product  (i.e.,  element  by  element  matrix  multipli¬ 
cation).  Then,  compute 


XAt  =  MFFTSat.  {XA}  = 


FFTS/Vpixai} 
FFTSaTp  {xaa} 

FFTS.  {xa^r} 


where 


XA  =  X  0  A  = 


and  FFTS^p  is  as  defined  before. 

r  ^ 

Step  3  Compute  the  iV-point  window  column-vector  g=  5r(l)  ^(2)  •••  g{N) 

and  create  G  (an  NxNp  matrix)  formed  by  Np  columns  of  g.  Also,  create  the 
exponential  matrix  E  as  follows: 


2 

••• 

2 

••• 

••• 
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where  = 

II 

1 

J 

•D* 

'  for  jfc  =  - 

At 

2  »  •  •  •  ’ 

2 

1  and  n  =  1 ,  •  •  • , 

X*  which  is 

defined  . 

as: 

■  +  1) 

+  1) 

:^'(f  +  l) 

X*  = 

x-(f  -f  2) 

+  2) 

■■■  x-{f  +  2) 

.  +  N) 

+  N) 

■■■  x*{^  +  N) 

Then,  create 


(30) 


where  x*(n)  is  the  complex  conjugate  of  a:(n).  Finally,  compute 


Xg  =XAt©E©X*, 


(31) 


also  an  NxNp  matrix. 


Step  4  The  cyclic  spectrum  in  matrix  form  (but  without  mapping)  is  found  to  be 


Sx  =  MFFTSiv  {Xg}  , 


(32) 


an  NxNp  matrix  where  each  element  is  denoted  for  —  l  and 


N„ 


-f  <9<f-l. 


Step  5  The  last  step  is  to  map  the  elements  of  the  matrix  Sx  onto  *?“(/),  i.e.,  Sx^,^  ^xif) 
where  a  and  /  are  related  to  k  and  ^  as  in  (24). 


As  can  be  seen,  this  is  a  highly  parallel  algorithm  suitable  for  implementation 
using  parallel  architectures.  However,  the  mapping  may  not  necessarily  be  a  trivial  task 
(see  next  section). 

Although  this  paper  does  not  intend  to  go  into  the  details  of  a  hardware  archi¬ 
tecture  using  VLSI  technology,  a  block  diagram  of  the  implementation  of  the  SSCA  is 
presented  in  Figure  1.  This  provides  a  succinct  and  complete  visual  form  of  the  SSCA. 


3.3  SSCA  Analysis 

The  analysis  of  the  SSCA  can  be  based  on  two  different  aspects:  the  computa¬ 
tional  requirements  and  the  accuracy  of  the  estimates.  The  problem  of  mapping  associated 
with  the  SSCA  is  treated  separately. 
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x(n):  block  of  N+N,  samples 


N/2  samples  N  s^ples  N/2  s|unpl 
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Table  1:  SSCA  Computational  Requirements. 


Operation 

#  of 

Cmplx  Mltply 

#of 

Cmplx  Addtn 

#  of 

Real  Mltply 

#of 

Real  Addtn 

Window  a{n) 

Np 

- 

2Np 

- 

Ap-point  FFT 

flog,  AT, 

Np  log2  Np 

2Np  log2  Np 

3  Ap  log2  Np 

Down-convert -fproduct 

2Np 

- 

8Np 

^Np 

Window  g{n) 

A 

- 

2N 

- 

A-point  FFT 

flog^A 

AlogjA 

2Alog2A 

3Alog2A 

3.3.1  Computational  Requirements 

Analyzing  the  SSCA  in  terms  of  computational  requirements  boils  down  to  an¬ 
alyzing  the  number  of  various  arithmetic  operations  required  for  each  steps  of  the  algo¬ 
rithm.  Note  however  that  the  analysis  is  made  with  respect  to  sequential  computations 
requirements. 

Assuming  a(r)  and  g{n)  are  already  defined  before  the  beginning  of  the  algo¬ 
rithm,  the  number  of  arithmetic  operations  required  for  each  step  (ignoring  the  mapping 
operation)  would  be  as  given  in  Table  1.®  Note  that  the  arithmetic  operations  under  in¬ 
vestigation  are  complex  multiplications  and  complex  additions  which  are  transposed  into 
real  multiplications  and  real  additions  in  the  two  last  columns. 

Denoting  the  total  number  of  real  multiplies  and  additions  for  the  single-processor 
as  and  respectively,  then  the  computational  requirements  are  found  to  be 

C:^^F{2  +  2\og,F)  (33) 

and 

C;.  =  f(4  +  31ogjf),  (34) 

where  F  =  N  •  NpJ  Now,  for  a  multi-processor  architecture  (although  not  necessarily 

®A  radix-2  FFT  algorithm  is  used. 

^For  the  single-processor,  the  first  three  operations  in  Table  1  are  performed  N  times  and  the  last 
two,  Np  times. 
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optimized),  the  computational  requirements  are 


CZ  =  N,{10  +  2  log2  N,)  +  N  (2 +  2  log,  N)  (35) 

and 

CZ  =  iVp  (4  +  3  log2  iVp)  +  iV  (3  log2  N)  (36) 

where  CZ  ^Jid  are  the  multi-processor  equivalent  of  CZ  ^'iid  (7*^.  Note  that  if  x{n)  is 
assumed  to  be  a  real  signal,  these  figures  will  go  down.  Clearly,  the  possibility  of  using  a 
parallel  architecture  with  the  SSCA  reduces  considerably  the  computational  requirements. 


3.3.2  Accuracy  of  the  Estimates 


The  accuracy  of  the  SSCA  for  estimating  the  cyclic  spectrum  is  related  to  the 
parameters  N,  Np,  and  the  type  of  data  tapering  window  chosen,  i.e.,  a(n)  and  g{n).  The 
frequency  resolution.  A/,  the  cycle  frequency  resolution,  Aq,  and  the  resolution  product, 
A/At  are  all  specified  in  terms  of  N  and  Np  as  stated  in  (18).  However,  the  choice 
of  a(n)  and  g{n)  will  impact  on  the  accuracy  of  the  point  estimates.  Essentially,  the 
bandwidth  of  a(n)  is  assumed  to  be  Aa  =  l/Np  (assuming  /,  =  1)  and  the  bandwidth  of 
g{n),  Ag  =  1/N.  The  kernel  transform  associated  with  the  SSCA  is  found  from  (17)  and 
is  given  by  [3]  [4]  [8] 

M(a,  /)  =  G(a  -  a„M  (/  -  /„  +  5^)  ,  (37) 


where  (?(•)  and  A(-)  are  the  Fourier  Transforms  of  g{n)  and  a(n)  respectively.  Also  recall 
that 


k  q 
f  _ l_ 

2Np  2N 


(38) 

(39) 


where  k  =  -Npf2, . . . ,  Np/2  —  1  and  q  =  —N/2, . . . ,  N/2  -  1  (again  assuming  fs  =  1).  As 
a  result,  A(^f  -  fo  +  is  a  strip  approximately  Aa  wide  along  the  line  a  =  2fk  —  2f 

(where  fk  =  k/Np).  This  implies  that  all  point  estimates  have  a  constant  frequency 
resolution  of  A/  =  Aa.  In  addition,  they  also  possess  a  constant  cycle  frequency  resolution 
which  is  determined  by  G{a),  i.e.,  Aa  =  1/N. 
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In  practice,  a(n)  or  g{n)  will  never  produce  perfect  bandpass  filters  and  some 
parameters  such  as  bandwidth  and  attenuation  are  always  cissociated  with  practical  filters. 
A  judicious  choice  for  the  best  windows  to  use  is  easier  to  make  by  experimentation.  The 
family  of  Kaiser  windows  allows  (via  the  parameter  P)  various  attenuations  and  therefore 
constitutes  a  good  choice  for  experimentation.  As  a  general  conclusion  coming  out  of  these 
investigations,  it  is  important  that  a{n)  has  a  large  attenuation  while  g{n)  has  a  small 
bandwidth  characteristic.  Since  both  attenuation  and  bandwidth  are  intimately  related, 
obtaining  a  larger  attenuation  is  equivalent  to  trading-off  the  bandwidth  (i.e.,  getting  a 
larger  bandwidth)  and  vice-versa.  As  a  result,  a{n)  will  have  a  large  attenuation  but  also 
a  large  bandwidth  and  g{n)  a  small  bandwidth  but  also  a  small  attenuation. 


3.3.3  Mapping  Associated  with  the  SSCA 

One  important  issue  of  the  SSCA  concerns  the  mapping  operation  as  briefly 
explained  in  Step  5.  This  will  now  be  discussed  in  more  depth. 

Before  going  into  the  details  of  the  algorithm  required  for  the  mapping,  the 
region  of  support  for  the  cyclic  spectrum  in  the  bifrequency  plane  is  first  examined.  This 
is  done  by  visualizing  the  mapping  of  k  and  q  onto  /  and  oc  in  Figure  2.  In  Figure  2(a), 
A  =  32  and  Np  =  S  are  used  and  therefore  k  =  —4,  . . . ,  3  and  q  =  —16,  . . . ,  15.  Each 
point  in  this  bifrequency  plane  represents  a  point  estimate  and  is  located  at  the  center  of 
a  rectangle  having  widths  Aa  =  1/Np  on  the  /-axis  and  Aa  =  1/N  on  the  a-axis.  Note 
that  /  is  minimum  at  9  =  15  (=  -y  —  1)  and  k  =  —4  (=  -^),  and  is  maximum  &t  q  =  —16 
(=  — y)  and  k  =  3  ^  —  1),  while  a  is  minimum  at  q  —  —16  and  k  =  —4,  and  is 

maximum  at  9  =  15  and  k  =  3.  As  can  be  seen,  each  negative-slope  diagonal  corresponds 
to  a  strip  for  a  particular  value  of  k  (the  bottom  diagonal  being  k  =  —  A/p/2  and  the  top 
one  k  =  Npf2  —  1).  Each  strip  contains  all  the  q  values  starting  from  the  top  of  any  strip 
at  q  =  —N/2  to  q  =  N/2  —  1  at  the  bottom.  Figure  2(b)  shows  the  mapping  for  higher 
N  and  A/p,  providing  smaller  Aa  and  Aa  (i.e.,  higher  frequency  resolutions).  Some  useful 
observations  are  noted:  (1)  all  the  point  estimates  having  the  same  value  of  a  will  form 
the  continuous  spectrum  related  to  this  particular  cycle  frequency  (recall  that  the  cyclic 
spectrum  is  discretely  distributed  along  the  cycle  frequency);  (2)  geometrically,  the  region 
of  support  is  found  to  be  a  diamond  (with  the  range  of  a  being  twice  that  of  /)  in  the 
bifrequency  plane;  and  (3)  this  diamond  is  not  perfect  in  the  sense  that  the  corners  of 
the  diamond  clearly  do  not  align.  The  coverage  on  the  negative  and  positive  frequency  is 
not  the  same  but,  nevertheless,  as  N  and  Np  goes  to  infinity,  it  tends  to  cover  almost  the 
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same  region  on  both  sides 


k  Md  q  MiB  f  nd  alpte  (N«32. 


Jyctal  fw^Macy.  f 


(a) 


MqipiBi  of  ktad  qflBiof  nd  i^la  (N*512.  Np*16) 


An  algorithm  for  implementing  the  mapping  operation  is  derived  in  Figure  3. 
To  start  with,  the  algorithm  assumes  that  Sx{q,k)  (from  Step  4)  has  been  computed 
and  is  available.  Then,  a  good  strategy  for  storing  that  data  for  any  further  use  (as  for 
plotting  or  analysis)  is  to  store  the  points  estimates  for  each  a.  Therefore,  the  output 
data  format  is  chosen  to  be  as  follows:  a  global  header,  followed  by  a  slice  header  (for 
each  a)  and  the  points  estimates  values  for  this  slice  with  this  latter  process  repeated 


14 


Global  Header 

surf_lype  (type  of  all  the  points  estimates  (1)  REAL  a  (2)  COMPLEX) 
niiin_slice/  Oob,  Omn  (#  of  slices,  the  minimum  and  maximum  value  of  a) 
max_nimi_f,  hm,  fun  (max  #  of  points  of  one  paiticulai  a, 

the  minimum  and  maximum  value  of  f  on  the  overall  gra;^) 

Slice#l 

a  (the  particular  value  of  a  for  this  slice) 

numj,  hmv  W  (#  of  pants  fa  tins  slice,  the  minimum  and  maximum  value 
off  for  this  slice) 

&((«■) ...  Sit(W)  (all  the  points  estimaies  fa  this  slice)  _ 


Slice  #num_slice 

a  (the  particular  value  of  a  fa  this  sice) 

numj,  hniv  fan  (#  of  points  fa  this  slice,  the  mininium  and  maximum  value 
tff  fa  this  slice) 

S»(fiiii) ...  S«(Uii)  (all  the  points  estimates  fa  this  slice) _ 

Figure  4:  Output  Data  File  Format. 
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for  all  slices  (see  Figure  4).®  The  algorithm  maps  the  points  estimates  using  this  format 
and  thus  enters  a  loop  for  each  value  of  a  (the  a-loop)  and  then  an  inner-loop  (the  f- 
loop)  for  each  value  of  /  possible  according  to  this  particular  a.  It  takes  into  account  all 
useful  geometrical  information  found  from  the  region  of  support  investigation  above  and 
represents  an  effective  way  of  mapping  the  data.  Note  however  that  if  visualization  of 
the  CS  is  not  required,  then  a  very  simple  skew  storage  mechanism  would  align  all  points 
estimates  with  respect  to  a  only.® 

3.4  Cyclic  Feature  Estimation  with  SSCA 

The  CS  constitutes  a  nice  tool  for  spectral  analysis  of  second-order  wide  sense 
cyclostationary  signals.  However,  it  is  often  necessary  to  only  identify  the  cyclic  features. 
In  fact,  all  signal  processing  methods  exploiting  cyclostationary  properties  of  signals  (eg 
cyclic  direction  finding  methods)  require  at  least  reliable  estimates  of  the  signal  of  inter¬ 
est’s  cyclic  features  in  order  for  them  to  be  useful.  One  can  derive  the  second-order  cyclic 
features  from  the  CS  estimate  such  as  the  one  produced  with  the  SSCA. 

With  the  SSCA,  the  number  of  point  estimates  in  the  CS  is  x  Np  but  the 
cycle  frequency  resolution  is  equal  to  ^  x  Hz.  The  maximum  number  of  points  used  for 
detecting  cyclic  features  should  therefore  be  no  greater  than  N  in  the  case  of  real  signals 
and  2N  —  ^  for  complex  signals.  Such  functions  are  defined  below. 

The  simplest  function  is  the  max-cut  function,  defined  by  the  author  to  be 

MAX(a)  ^  max  [|5"(/)P]  .  (40) 

It  retains  the  maximum  magnitude  squared  of  each  a-cut.  The  max-cut  function  can 
easily  be  modified  into 

which  is  hereinafter  referred  to  as  the  normalized  log  max-cut  (NLMC)  function.  It  is 
easily  seen  that 

LM(0)  =  0  (42) 

®This  output  data  file  format  is  used  as  an  input  by  the  plot^xaf  program  to  plot  the  cyclic  spectrum. 
This  software  is  sold  by  Statistical  Signal  Processing  Inc.  (SSPI). 

®This  could  be  the  case  for  example  when  hardware  implementation  is  required  and  that  the  cyclic 
spectrum  is  used  to  detect  cyclic  features. 
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and  furthermore 

LM(a)  <  0  (43) 

since  MAX(O)  should  be  greater  to  or  equal  to  MAX(q;).  Even  though  these  functions 
could  be  successfully  used  in  practice  to  detect  cyclic  features,  the  best  theoretical  function 
(or  measure)  is  called  the  cycle-frequency  decomposed  measure  of  degree  of  cyclostation- 
arity  and  is  given  by  [10] 


IZomr)fdT 

irLmr)fdT- 


(44) 


An  estimate  of  this  measure  found  from  the  CS  estimate  given  by  the  SSCA  is  defined 
here  to  be 

>  A  ,  ,  Ey|S?(/)l’ 

Dia)  =  wla)  - — - — y  {45) 

E/IW)I 


where  w(a)  denotes  a  weighting  function  of  a.  This  weighting  function  takes  into  consid¬ 
eration  the  reliability  of  the  CS  measurement  at  a.  Using 


numf(a) 

numf(0) 


(46) 


with  numf(a)  denoting  the  number  of  point  estimates  computed  at  a,  the  effect  of  the 
varying  number  of  points  along  a  introduces  a  false  relationship  between  the  cyclic  feature 
strengths.  Specifically,  going  from  a  =  0  to  a  =  ±1,  the  cyclic  features  strength  would 
decrease  even  if  they  should  be  equal  in  reality.  A  way  of  avoiding  this  is  to  make  sure  the 
signal  does  not  have  any  frequency  component  higher  than  /s/4  or  lower  than  —/s/4  and 
consider  the  CS  cuts  only  from  a  =  —fs/2  to  a  =  /a/2,  therefore  making  the  number  of 
points  to  be  added  for  each  cut  constant.  This  is  equivalent  to  interpolating  the  signal  by 
two  and  considering  only  the  SSCA  points  estimates  lying  within  the  rectangular  region 
/a/2  <  a  <  fa/2,  —/a/4  <  /  <  /a/4.  Consequently,  u;(q:)  =  1.  This  method  however  will 
double  the  number  of  points  to  be  considered  in  the  SSCA  and  complicates  the  realization 
of  the  algorithm.  The  author  has  tried  another  way  of  estimating  DCS“  with  the  SSCA 
points  estimates.  Since  DCS"  is  resulting  from  a  ratio  of  infinite  integrals,  it  could  easily 
be  seen  as  a  ratio  of  the  averaged  power  at  each  cut  and  the  averaged  power  at  a  =  0. 
To  realize  this,  we  simply  have  to  set 


numf(0) 

numf(Q) 


(47) 
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Nevertheless,  this  estimate  suffers  from  reliability  problems  as  we  approach  a  =  ±/s  since 
less  points  are  taken  into  account  in  the  averaging  operation  and  that  the  feature  will  get 
amplified  by  noisy  points.  As  a  result,  the  estimate  using 

w(a)  =  1  (48) 

seems  to  be  the  compromise.  From  a  cyclic  feature  detection  point  of  view,  this  last 
estimation  method  will  be  adequate  and  reliable  but  if  someone  is  looking  for  feature 
strength  relationships,  a  method  similar  to  the  one  described  above  and  using  interpolation 
may  be  required. 

Again,  this  measure  can  be  represented  in  dB  by  taking  the  logarithm,  i.e. 

LD(a)  =  101ogioD(a)  (49) 

or  equivalently 

LD(a)  =  101og,o  El^x(/)I'  -C-  (50) 

.  / 

where  C  is  given  by 

c=ioiog.„  Ei-^(/)r 

.  / 

This  measure  provides  a  more  reliable  function  for  studying  the  signals  cyclic  features 
distribution  and  should,  as  the  NLMC  function,  always  be  smaller  or  equal  to  0  (LD(q;)  < 
0).  As  a  result,  LD(a),  hereinafter  referred  to  as  the  Cyclic  Feature  Function  (CFF),  is 
preferred. 

4.0  RESULTS  OF  THE  SSCA  ENCODED  IN  C 


The  algorithm  explained  in  the  last  section  has  been  encoded  in  the  C  language 
on  a  SUN  IPX  computer.  The  program  is  called  ssca  and  its  usage  is  given  as  follows: 

ssca  [options]  filenamel  [filename2] 

The  filenamel  is  the  name  of  the  source  file  containing  the  block  of  samples  of 
the  signal  to  be  analyzed.  The  second  optional  file  name,  filename2,  is  required  only  when 
the  cross-correlation  cyclic  spectrum  is  computed.  The  different  options  that  can  be  used 
are  listed  in  Table  2. 
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Table  2:  Options  for  the  ssca  Program. 


Options 

Description 

-a  value 

specifies  the  amount  of  attenuation  (in  dB)  for  the  1st  Kaiser  window  (Np 
samples  long)  -  [Default  96.0  dB] 

-b  value 

specifies  the  amount  of  attenuation  (in  dB)  for  the  2nd  Kaiser  window  (N 
samples  long)  -  [Default  =>  20.0  dB] 

-C 

computes  the  conjuguate  cyclic  spectrum 

-c 

computes  the  cross-correlation  cyclic  spectrum  (filename2  must  be  specified  in 
this  case) 

-f  char 

specifies  the  format  of  the  data  to  be  produced  in  the  output  files  (if  char  =  a, 
the  output  is  ASCII  and  if  char  =  b,  it  is  binary)  -  [Default  b] 

-help 

displays  information  on  how  to  obtain  on  line  help 

-M 

in  addition  to  the  cyclic  spectrum  output,  also  outputs  the  maximum  value  of 
the  cyclic  spectrum  for  each  cycle  frequency  in  a  separate  output  file  (using 
extension  .max) 

-N  value 

specifies  the  number  of  data  samples  being  analyzed  (i.e.,  N  samples)  -  [Default 
1024] 

-P  value 

specifies  the  number  of  data  samples  used  in  the  channelizer  (i.e.,  Np  samples) 

-  [Default  8] 

-o  file 

specifies  the  output  file  name  -  [Default  surf_out] 

-d 

in  addition  to  the  cyclic  spectrum  output,  also  outputs  the  degree  of  cyclo- 
stationarity  for  each  cycle  frequency  in  a  separate  output  file  (using  extension 
.int) 

-z 

in  addition  to  the  cyclic  spectrum  output,  also  outputs  the  value  of  the  cyclic 
spectrum  at  the  zero  frequency  (or  closest  to  zero)  for  each  cycle  frequency  in 
a  separate  output  file  (using  extension  .zero) 

-s  value 

decimates  the  number  of  output  spectrum  cut  to  value  cuts  without  affecting 
the  resolution.  -  [Default  =>  1024] 
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Table  3:  Signals  Parameters  (for  fs  =  lOOkHz). 


Signal 

Keying  Rate 

Carrier  Frequency 

Pulse  Type  q{t) 

Bandwidth 

Power  (dBW) 

B 

4800  bauds 

23  kHz 

rect  BL 

4*4800 

0.0 

B 

10000  bauds 

23.5  kHz 

rect  BL 

3*10000 

0.0 

C 

9000  bauds 

29.5  kHz 

rect  BL 

2*9000 

0.0 

NO 

— 

— 

— 

fsl2 

0.0 

N5 

— 

— 

— 

fs/2 

5.0 

Using  the  ssca  program,  this  section  will  compare  cyclic  spectrum  results  with 
theoretical  expectations.  Note  that  the  program  used  to  display  the  3-dimensional  spec¬ 
trum  is  called  plot_sxaf  and  is  sold  by  Statistical  Signal  Processing  Inc  along  with  another 
type  of  cyclic  spectral  estimation  algorithm  using  direct  frequency-smoothing. 


4.1  Signal  Descriptions 

The  type  of  signal  used  for  the  analysis  is  the  popular  Phase-Shift  Keying  (PSK) 
modulated  signal.  Specifically,  binary  PSK  (BPSK)  modulated  signals  are  used  with 
additive  white  Gaussian  noise  (AWGN).  A  BPSK  signal  can  be  mathematically  described 
in  terms  of  pulse-amplitude  modulation  (PAM)  as 

z(t)  =  a{t)  cos  (27r/ct  -|-  (/>o)  ,  (52) 

where  the  signal  envelope  a{t)  is  given  by 

00 

o(0  =  2  -  TiTk  -  to)  ■  (53) 

n=— oo 

Note  that  a„  is  a  binary  sequence  (±1),  q{t)  is  a  finite  energy  keying  envelope,  and 
Tk  =  I/Ik  is  the  symbol  period  expressed  in  seconds/symbols.^®  Table  3  and  Table  4 
describe  three  BPSK  and  two  noise  signals  under  investigation,  along  with  their  specific 
parameters. 

^°Note  that  1  baud  corresponds  to  1  symbol/sec. 
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Table  4:  Normalized  Signals  Parameters  (/,  =  IHz). 


4.2  Theoretical  Expected  Results 


The  theoretical  CS  of  a  BPSK  signal  is  well  known  [11]  and  is  derived  below. 


We  have 


R“(r)  =  (a;(t  +  r/2)x(i  -  T/2)e-^2’^“‘) 

=  (a{t  +  t/2)  cos  [2Tr fc{t  +  t/2)  +  ^o]  a(<  —  t/2)  cos  [2iTfc{t  —  r/2)  +  <f>Q] 

=  (a{t  +  Tf2)a{t  -  'r/2)i  [cos{2ir Jct)  +  cos(47r/c<  +  2<f>o)] 

=  +  r/2)a(t  -  r/2)  cos(27r/cr)e“''^"*^ 

+  +  T/2)a(t  -  r/2)  cos(4^/cf  +  2<^o)e"-'^“‘^ 


/  1  ^i(4ir/ct+2^o)  _L  p-j(47r/c«+2.^o)  \ 

+  ^ia(t  +  r/2)a(t-r/2).^ - ± - 

=  •  cos(27r/cr) 

+i  ^a(i  +  r/2)a(<  -  r/2)e-^='’^<“-2^'>‘e^’2^'>  +  a(t  +  r/2)a(t  -  r/2)e-^'2^(“+2^=)‘e-^‘2^‘>) 
=  •  cos(27r/,r)  +  (54) 

where  i?“(r)  is  the  cyclic  autocorrelation  of  the  PAM  signal  a{t).  Using  the  cyclic  Wiener 
relation,  one  can  find  the  CS  of  x{t)  by  Fourier  transforming  R^{t),  i.e., 


s^if)  =  ^{r:{t)} 
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=  \  [•?:(/  -  fc) + s:if + fc) + sr'^Hf)  •  . 


For  a{t)  as  expressed  in  (53),  its  CS  is  found  to  be 


iQ(/  +  a/2m/-a/2)S,(/  +  Q/2)e--’=’“‘»  ,  a  =  p/rK 
(/)  =  {  W°) 

I  0  ,  otherwise 

where  p  is  an  integer,  Q{f)  =  ^  {^(t)},  Saif)  =  5Zjb=-oo  Raik)  exp{-j2TrkfTK)  and  Raik)  = 
<  ajt+„a„  >  (discrete-time  average  of  the  sequence  over  n).^^  Using  these  results,  we  can 
then  express  S°{f)  as 


i[W-/2W(/-/i)s.(/-W 

+<3(/  +  /.)«•(/  +  h)s.(f  +  /.)]  ■ 


,  a  =  p/Tk 


stif)  = 


[Qif  -  h)QV  +  h)Uf  - 

+^(/  +  h)QV  -  h)Saif  +  ,  a  =  ±2/,  pITk 

(57) 


where 


/i  —  /c  +  q:/2 
/2  =  /c-a/2. 

This  is  the  general  expression  for  the  CS  of  a  BPSK  modulated  signal  as  expressed  in  (52) 
and  (53). 

The  actual  choice  of  g(t)  (i.e.,  a  band-limited  rectangular  pulse)  is  given  by 


or  equivalently 


qit)  =  rectr^(t)  ♦  ^(t) 


<?(/)  =  ■  HUf) 


where  rect7’^(<)  is  defined  as  being  the  rectangular  pulse  of  time  duration  Tk  and  hip{t)  is 


^^Note  that  for  {cn}  =  ±1  and  independently  identically  distributed  (i.i.d.),  Saif)  =  1- 
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the  impulse  response  of  the  filter  band-limiting  the  signal  having  Hip{f)  as  its  frequency 
response.  Replacing  Q{f)  in  the  expression  for  the  CS  gives 


g— >2>rato 

*Tk 


smMf-h)TKWnMf-ft) 


■UM  -  h)HM  -  fi)SM  -  M 

I  8m(7r(/+/i)ry)-sin(^r(/+/2)rff) 

'  5r(/+/i)-ir(/+/2) 

■HM  +  MHUf  +  h)S.U  +  /i)] 


,  a  =  pITk 


cat  f\  _  ,  e-i2Tato  \am(Tr(f-j2)TK)-sm(Tr(f+f2)Tjc) 

^ ATk  [  Af-j2)-AJ+h) 

■HUf  -  h]HUf  +  f2)S,(f  -  h)  ■ 

■  «in(^r(/+/i)Ty)-gin(7r(/-/i)ry) 

'  ^U+h)-Af-h) 

■Hlpif  +  h)HM  -  h)S.(f  +  /i)  •  ,  a  =  ±2/.  +  p/Tk 


,  otherwise 


and  furthermore,  choosing  the  sequence  {a(n)}  of  ±1  to  be  independently  identically 
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distributed  (i.i.d.)  causes  Sa(f)  =  1,  and  therefore  5“(/)  reduces  to 


g->2jrQ«o 

^Tk 


8in(^r(/-/2)rK)-Mn(^(/-/i) 


-  h)Hi,{f  -  h) 

I  «ip(’r(/+/i)rff)-8in(ir(/+/2)ry) 
’r(/+/i).jr(/+/2) 

^  h)Hi,{f  +  h)] 


e-j2irato  r8in(8-(/-/2)ry)-8in(ir(/+/2)ry) 

4Tk  [  Af-h)-<J+h) 

+  h)  ■  ei(^-/^‘o+2^o) 

I  MP(^(/4-/i)rAr)-Bin(7r(/-/i)rK-) 


,  a  =  p/Tii: 

,  a  =  ±2/c  +  p/Tx 

,  otherwise 


(61) 


This  last  expression  shows  that  the  BPSK  will  have  spectral  correlation  components  at 
a  =  p/Tk  and  a  =  ±2/c  +  pfT^  (for  any  integer  p).  In  addition,  the  CS  should  look 
like  a  sine-squared  function  at  a  =  0  (/  =  ±fc)  and  a  =  ±2/c  (/  =  0),  band-limited  by 
the  bandwidth  of  the  filter.  For  p  ^  0,  the  CS  more  or  less  corresponds  to  spectral 
correlation  due  to  the  p**  lobe  of  the  spectrum.  Since  the  number  of  lobes  is  band-limited 
to  the  bandwidth  of  the  filter,  the  number  of  spectral  correlation  terms  will  be  reduced 
accordingly.  The  expected  cyclic  features  for  each  signal  as  described  in  Table  4  are  shown 
in  Table  5.^^ 


4.3  Experimental  Results 

In  this  section,  the  CSs  of  the  various  generated  signals  (A,  B,  C,  NO,  and  N5) 
are  computed  using  the  ssca  program.  The  effects  of  varying  some  parameters  in  the 
algorithm  are  investigated.  All  CS  3D- plots  only  show  |5'"(/)|  on  a  linear  scale.  The  /- 
axis  and  a-axis  appear  on  the  left  and  right  side  of  the  graph  respectively.  Also  note  that, 

the  table,  /  denotes  the  center  frequency  at  which  a  cyclic  feature  is  present. 
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Table  5:  Theoretical  Signals  Cyclic  Features  Distribution. 


Cyclic 

Feature 

A 

B 

C 

a  =  0 

/  =  ±0.23 

/  =  ±0.235 

/  =  ±0.295 

a  =  ±pITk 

a  =  ±0.048  •  p 

f  =  ±0.23 

a  =  ±0.1  •  p 

f  =  ±0.235 

a  =  ±0.09  •  p 

f  =  ±0.295 

a  =  ±2fc 

a  =  ±0.46 

/  =  o 

a  =  ±0.47 

/  =  o 

a  =  ±0.59 

/  =  o 

a  =  ±2/c  -1-  pITk 

a  =  ±0.46  ±  0.048  •  p 

/  =  o 

a  =  ±0.47  ±0.1  -p 

/  =  o 

a  =  ±0.59  ±  0.09 -p 

/  =  o 

since  the  option  -s  512  has  been  used  to  reduce  the  large  number  of  points  generated 
by  the  algorithm,  the  number  of  output  cuts  has  been  decimated  to  512  cuts  for  each 
reported  result. 

The  effect  of  varying  the  observation  time  (N)  on  the  CS  of  signal  A  is  shown 
in  Figure  5.  First,  increasing  the  observation  time  increases  the  strength  of  the  spectral 
correlation  terms,  provided  the  signal  is  present  for  the  entire  observation  time  (which 
is  the  case  here).  Second,  since  Np  is  kept  fixed  in  this  case,  increasing  N  corresponds 
to  increasing  At  A/  (=  N/Np).  Figure  5(a)  clearly  shows  undesired  spectral  components 
(due  to  measurement  unreliability).  Starting  from  (c)  (AtA/  =  16384/64  =  256),  the  CS 
appears  free  of  these  unwanted  components  and  therefore  seems  to  be  a  reliable  measure¬ 
ment  (qualitatively). 

Even  though  the  CS  3D-plot  is  a  useful  for  spectral  analysis  of  cyclostationary 
signals,  another  plot  derived  from  the  CS  and  more  suitable  for  cyclic  features  detection 
has  been  defined  in  an  earlier  section.  This  function,  the  CFF  (=LD(q;)),  estimates  the 
cycle  frequency-decomposed  measure  of  the  DCS  or  equivalently  represents  the  cyclic 
feature  strengths  in  dB  distributed  over  a.  This  function  was  produced  with  the  ssca 
program  by  using  the  option  -d.  As  shown  in  Figure  6,  the  CFF  graphs  (using  the  same 
parameters  conditions  as  the  previous  figure)  provide  better  pictures  for  studying  the 
cyclic  features  along  The  spectral  shape  however  is  lost.  In  this  section,  the  CFF 

^^The  CFF  graphs  show  a  dashed  line  arbitrarily  set  to  -40  dB  representing  the  minimum  value  for  a 


26 


-10 


-10 


(a)  A  (A=1024,Ap=64)  (b)  A  (A=4096,  Arp=64) 


(c)  A  (A= 16384,  JVp=64) 


(d)  A  (iV=32768,  iVp=64) 


Figure  6:  CFF  of  a  BPSK  Signal:  Varying  the  Observation  Time. 


28 


is  used  in  order  to  more  easily  appreciate  the  distribution  of  cyclic  spectral  correlation 
components  along  the  cycle  frequencies.  Each  feature  reveals  itself  as  a  discrete  peak  at 
a  specific  a  corresponding  to  the  theoretical  value  mentioned  earlier  for  signal  A.  The 
measurement  noise  is  visible  on  each  graph  and  clearly  decreases  from  (a)  to  (d)  such 
that  peaks  exhibit  higher  SNR. 

Next,  the  effect  of  varying  the  frequency  resolution  A/  =  1/Np  on  signal  A  is 
studied  in  Figure  7  and  Figure  8.  Using  A/  =  and  |  in  (a),  (b),  (c)  and  (d) 

respectively,  a  loss  in  frequency  resolution  appears  in  the  CS  (from  (a)  to  (d)).  On  the 
other  hand,  the  product  AtAf  is  increased.  However,  there  should  be  some  constraints 
imposed  on  A/,  i.e.,  for  a  BPSK  signal  having  keying  rate  1/Ta, 

A/ <5^.  (62) 

In  Figure  8(d),  for  example.  A/  =  |  =  0.125  and  1/2T^  =  0.048/2  =  0.024  and  therefore 
Af  >  1/2Ta.  As  can  be  seen,  this  situation  turns  out  to  be  undesirable  since  the  baud  rate 
feature  decreases  in  strength  and  a  highly  distorted  relationship  between  cyclic  features 
arises. 


Figure  9  to  Figure  12  demonstrate  the  same  varying  parameters  as  for  signal  A 
but  applied  to  a  WGN  signal  NO.  As  observed,  increasing  the  observation  time  without 
changing  the  frequency  resolution  reduces  the  noise  measurement  and  its  variance.  As 
seen,  the  number  of  false  detected  features  greatly  reduces  from  (a)  to  (d).  In  general,  the 
CFF  inherently  takes  a  specific  shape  such  that  the  function  clearly  depicts  a  higher  noise 
measurement  variance  as  we  approach  a  =  ±1  and  the  feature  strengths  decrease  from 
(a)  to  (d).  This  particular  shaping  of  the  CFF  is  due  to  the  region  of  support  which  was 
shown  to  be  a  diamond  and  therefore  contains  fewer  points  to  rely  on  at  the  extremes. 
The  effect  of  reducing  the  frequency  resolution  is  studied  in  Figure  12.  Due  to  the  signals 
under  investigations,  N  and  Np  have  been  chosen  to  be  32768  and  64  respectively  for  all 
subsequent  figures. 

The  CS  and  CFF  of  all  noise-free  BPSK  signals  are  shown  in  Figure  13  and 
Figure  14  respectively.  In  (d),  all  signals  have  been  combined  and  despite  the  fact  that 
it  is  not  possible  to  detect  the  3  signals  and  their  parameters  by  only  analyzing  the 
conventional  PSD  (at  a  =  0),  investigations  of  the  CS  or  CFF  reveal  every  signals’ 

peak  to  be  considered.  The  dotted  lines  on  the  graphs  appear  on  peaks  which  have  been  detected  using 
a  detection  algorithm  derived  by  the  author.  Quantitative  results  are  provided  in  Appendix  A. 
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(c)  A  (Ar=32768,  Np=16) 


Figure  7:  CS  of  a  BPSK  Signal: 


(b)  A  (A=32768,  Np=Z2) 


(d)  A  (Ar=32768,  Np=S) 


arying  the  Frequency  Resolution. 


(a)  NO  (Ar=1024,  Np=64) 


(b)  NO  (Ar=4096,  Np=64) 


Figure  9:  CS  of  Noise:  Varying  the  Observation  Time. 


(a)  NO  (iV=1024,  Np=M) 


(c)  NO  (Ar=16384,  JVp=64) 


Figure  10:  CFF  of  Noise: 


FMture  Strength  (dB) 


(a)  NO  (7V=32768,  Np=64) 


(c)  NO  (iV=32768,  Np=16) 


Figure  11:  CS  of  Noise: 


(b)  NO  (Ar=32768,  Np=32) 


(d)  NO  (iV=32768,  Np=8) 


the  Frequency  Resolution. 


Feature  Strength  (dB)  Feature  Strength  (dB) 


(a)  NO  (iV=32768,  JVp=64)  (b)  NO  (Ar=32768,  Np=i2) 


(c)  NO  (iV=32768,  Np=16) 


(d)  NO  (Ar=32768,  Np=8) 


Figure  12:  CFF  of  Noise:  Varying  the  Frequency  Resolution. 
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parameters  (see  Table  10  of  Appendix  A).  The  effect  of  adding  noise  to  these  signals  is 
studied  from  Figure  15  to  Figure  18.  Again,  as  the  results  tabulated  in  Appendix  A 
quantify,  all  major  features  can  still  be  detected.  This  undoubtly  demonstrates  the  high 
degree  of  robustness  and  precision  of  the  CS  (or  CFF)  tools  for  analyzing  cyclostationary 
communications  signals. 

4.4  SSCA  Benchmarks  on  a  SUN  IPX 

The  computation  time  required  for  a  SUN  IPX  to  perform  the  SSCA  using 
specific  parameters  are  reported  in  Table  6  for  reference.  All  timings  are  accurate  to 
within  1  second  and  do  not  include  the  first  and  last  steps  of  the  algorithm  which  read 
the  data,  form  the  windows  coefficients,  and  map  and  store  the  output.  Note  however 
that  the  time  needed  for  these  steps  is  relatively  negligible. 

Table  6:  Results  of  SSCA  Benchmarks  on  a  SUN  IPX. 


N 

N, 

AtAf 

Time  (sec) 

1024 

8 

128 

<  1 

1024 

16 

64 

1 

1024 

32 

32 

2 

1024 

64 

16 

4 

4096 

64 

64 

14 

8192 

64 

128 

30 

16384 

64 

256 

77 

32768 

64 

512 

169 

32768 

32 

1024 

81 

32768 

16 

2048 

41 

32768 

8 

4096 

20 

38 


(a)  A+NO 


(c)  C+NO 


Figure  15;  CS  of  Various  BPSK  Si 


(d)  A+B+C+NO 


s  in  Noise  NO  Using  iV=32768,  iVp=64. 


•1  ^.8  -0.6  -0.4  -0.2  0  -1  -0.8  -0.6  -0.4  -0.2  0 

Cycle  Frequency  (Hz)  Cyde  Frequency  (Hz) 


(a)  A+NO  (b)  B+NO 


-1  -0.8  -0.6  -0.4  -0.2  0  -1  -0.8  -0.6  -0.4  -0.2  0 

Cycle  Frequency  (Hz)  Cycle  Frequency  (Hz) 


(c)  C-fNO  (d)  A-fB-fC+NO 

Figure  16:  CFF  of  Various  BPSK  Signals  in  Noise  NO  Using  iV=32768,  7Vp=64. 
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(a)  A+N5 


(b)  B+N5 


(c)  C+N5 


(d)  A+B+C+N5 


Figure  18:  CFF  of  Various  BPSK  Signals  in  Noise  N5  Using  iV=32768,  iVp=64. 
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5.0  SUMMARY  AND  CONCLUDING  REMARKS 


This  report  discussed  the  implementation  of  an  efficient  algorithm,  the  SSCA,  for 
estimating  the  entire  cyclic  spectrum  of  discrete-time  domain  signals.  Following  a  review 
on  the  theory  of  spectral  correlation  for  cyclic  spectral  analysis  of  discrete-time  signals, 
the  SSCA  was  derived.  The  detailed  description  of  the  SSCA  showed  that  it  mainly 
reduces  in  complexity  to  a  two-dimensional  FFT  of  size  NxNp  and  that  it  is  suitable  for 
parallel  architectures.  Parameters  N  and  Np  specify  a  number  of  CS  estimate  attributes 
such  as  the  frequency  resolution  A/,  the  cycle  frequency  resolution  Aa,  and  the  time- 
frequency  resolution  product  A/At.  All  the  steps  required  to  implement  the  SSCA  were 
clearly  explained.  An  analysis  of  the  SSCA  revealed  its  computational  requirements.  The 
accuracy  of  the  estimates  were  also  considered  and  an  algorithm  for  mapping  the  SSCA 
output  was  provided.  Additionally,  the  cyclic  feature  function  (CFF)  was  defined  as  a 
means  to  detect  the  cyclic  features  from  the  SSCA. 

Results  of  the  SSCA  encoded  in  C  were  then  reported.  Theoretical  CSs  were 
derived  for  three  BPSK  signals  using  various  parameters  and  two  additive  white  Gaus¬ 
sian  noise  (AWGN)  signals.  Plots  of  the  CS  and  CFF  of  various  combinations  of  these 
signals  along  with  various  algorithmic  parameters  were  realized  with  the  ssca  program 
and  presented  to  the  reader.  Specifically,  the  effect  of  varying  the  observation  time  and 
the  frequency  resolution  were  observed  on  one  of  the  BPSK  signal  and  one  of  the  AWGN 
signal  independently.  Then,  the  CS  and  CFF  of  combinations  of  the  signals  under  inves¬ 
tigation  were  studied  while  all  SSCA  parameters  were  kept  fixed.  Cyclic  features  were 
detected  based  on  the  CFF.  These  experimental  results  show  that  all  major  cyclic  fea¬ 
tures  corresponding  directly  or  indirectly  to  signal  features  (e.g.  carrier  frequency  and 
baud  rate  in  this  case),  were  successfully  detected  even  in  highly  corrupted  conditions 
where  conventional  analysis  would  have  failed  to  discriminate  the  signals.  Finally,  some 
benchmarks  on  the  SUN  for  the  SSCA  were  provided  for  reference. 

In  conclusion,  the  CS  and  CFF  estimated  with  the  SSCA  have  proven  to  be 
valuable  tools  for  analyzing  second-order  cyclostationary  communication  signals  and,  by 
making  extensive  use  of  the  FFT,  to  provide  reliable  and  accurate  results  more  efficiently 
than  typical  CS  direct  estimation  methods.  Based  on  these  results,  it  is  anticipated  that 
the  SSCA  will  be  coded  in  a  multi-C40  processors  board  in  a  near  future  and  will  be 
processing  off-air  data.  Cyclic  spectral  analysis  of  live  data  composed  of  various  modu¬ 
lation  types  signals  will  then  be  performed  and  the  practical  usefulness  of  the  SSCA  will 
be  assessed. 
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A.O  CYCLIC  FEATURE  DETECTION  RESULTS 

These  results  were  found  based  on  the  Cyclic  Feature  Function  (CFF)  as  defined 
in  the  main  body  of  the  report.  A  cyclic  feature  detection  algorithm  has  been  applied  to 
the  CFF  producing  tables  of  detected  cyclic  features.  Each  table  shows  the  cyclic  feature 
value  a,  its  strength  in  dB,  and  to  which  signal  feature  it  corresponds.  The  term  FALSE 
means  that  the  feature  does  not  correlate  with  any  signal’s  cyclic  features.  The  detection 
algorithm  used  is  not  described  here  but  is  a  fairly  straightforward  manner  of  detecting 
the  features  based  on  the  CFF.  Only  results  of  the  SSCA  using  N  =  32768  and  Np  =  64 
are  tabulated. 


Table  7:  A  (Y=32768,  Aip=64),  21  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.04800 

-11.31 

-0.04800 

-11.31 

fkA 

0.09601 

-18.17 

-0.09601 

-18.17 

Vk, 

0.14401 

-33.48 

-0.14401 

-33.48 

3A. 

0.31601 

-35.39 

-0.31601 

-35.39 

2/c^  ~  3/a;^ 

0.36401 

-22.83 

-0.36401 

-22.83 

2/c/l  2/fc^ 

0.41199 

-16.43 

-0.41199 

-16.43 

2/ca  “  fkA 

0.45999 

-4.08 

-0.45999 

-4.08 

^fcA 

0.50800 

-14.34 

-0.50800 

-14.34 

^fcA  d"  fkA 

0.55600 

-20.17 

-0.55600 

-20.17 

^fcA  d" 

0.60400 

-34.42 

-0.60400 

-34.42 

2/cy|  d-  3/fc^ 
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Table  8:  B  (7V=32768,  Np=6i),  15  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.10001 

-10.00 

-0.10001 

-10.00 

fkg 

0.20001 

-25.80 

-0.20001 

-25.80 

2A, 

0.26999 

-28.36 

-0.26999 

-28.36 

2/cb  -  3/itB 

0.37000 

-12.81 

-0.37000 

-12.81 

2/cb  -  fkB 

0.47000 

-3.03 

-0.47000 

-3.03 

2/cb 

0.57001 

-13.26 

-0.57001 

-13.26 

2/cb  d*  fkB 

0.67001 

-29.30 

-0.67001 

-29.30 

2/cB+2/fcB 

Table  9:  C  {N—3276S,  Np=64:),  9  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.09000 

-10.24 

-0.09000 

-10.24 

fkc 

0.50000 

-13.06 

-0.50000 

-13.06 

2/cc  -  fkc 

0.59000 

-3.20 

-0.59000 

-3.20 

Vcc 

0.67999 

-13.84 

-0.67999 

-13.84 

^fcc  d"  fkc 

A-2 
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Table  10:  A+B+C  (iV=32768,  Ap=64),  35  Features  Detected. 


a 

Value  (dB) 

Q 

Value  (dB) 

0.00000 

0.00 

0.04800 

-18.15 

-0.04800 

-18.15 

0.09000 

-18.61 

-0.09000 

-18.61 

0.09601 

-24.83 

-0.09601 

-24.83 

0.10001 

-19.25 

-0.10001 

-19.25 

0.20001 

-35.28 

-0.20001 

-35.28 

0.27000 

-37.50 

-0.27000 

-37.50 

0.36401 

-29.65 

-0.36401 

-29.65 

0.37000 

-21.84 

-0.37000 

-21.84 

0.41199 

-23.10 

-0.41199 

-23.10 

0.45999 

-10.75 

-0.45999 

-10.75 

0.47000 

-12.26 

-0.47000 

-12.26 

0.50000 

-20.49 

-0.50000 

-20.49 

0.50800 

-21.21 

-0.50800 

-21.21 

0.57001 

-22.37 

-0.57001 

-22.37 

0.59000 

-11.61 

-0.59000 

-11.61 

0.67001 

-39.48 

-0.67001 

-39.48 

0.67999 

-22.36 

-0.67999 

-22.36 

Feature 

0 

fkA 

fkc 

VkA 

fkg 


2/fcB 

2/c. 

-Vks 

2/cyl 

-VkA 

2/cb 

-  fkB 

2/cx 

~  fkA 

^fcA 

Vcs 

2/cc 

-  Ike 

2/c^ 

+  fkA 

2/cb 

+  fkB 

2/cc 

2/cb 

+  2/*;b 

^fcc 

+  fkc 
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Table  11:  A+NO  (A=32768,  Np=M),  15  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.04800 

-13.53 

-0.04800 

-13.53 

f^A 

0.09601 

-21.15 

-0.09601 

-21.15 

VkA 

0.36401 

-25.66 

-0.36401 

-25.66 

VcA  -  VkA 

0.41199 

-18.59 

-0.41199 

-18.59 

^fcA  ~  f^A 

0.45999 

-6.57 

-0.45999 

-6.57 

^fcA 

0.50800 

-16.66 

-0.50800 

-16.66 

^fcA  d"  fkA 

0.55600 

-22.35 

-0.55600 

-22.35 

VCA  +  VkA 

Table  12:  B+NO  (A=32768,  Ap=64),  14  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.10001 

-13.57 

-0.10001 

-13.57 

fkg 

0.20001 

-27.66 

-0.20001 

-27.66 

Vks 

0.37000 

-16.05 

-0.37000 

-16.05 

^fcB  ~  fkB 

0.47000 

-6.69 

-0.47000 

-6.69 

VcB 

0.57001 

-17.10 

-0.57001 

-17.10 

^fcB  d"  fkB 

0.66999 

-30.82 

-0.66999 

-30.82 

FALSE 

0.73911 

-34.07 

FALSE 

A-4 
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Table  13:  C+NO  (A^=32768,  Np=6A),  9  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.09000 

-13.33 

-0.09000 

-13.33 

fkc 

0.50000 

-15.96 

-0.50000 

-15.96 

^fcc  ~~  fkc 

0.59000 

-6.27 

-0.59000 

-6.27 

2/cc 

0.67999 

-16.90 

-0.67999 

-16.90 

2/cc  +  fkc 
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Table  14:  A+B+C+NO  (^=32768,  Ap=64),  31  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.04800 

-19.28 

-0.04800 

-19.28 

f^A 

0.09000 

-19.35 

-0.09000 

-19.35 

fkc 

0.09601 

-25.77 

-0.09601 

-25.77 

VkA 

0.10001 

-20.54 

-0.10001 

-20.54 

fkB 

0.36401 

-30.09 

-0.36401 

-30.09 

^fcA  ~  "^fkA 

0.37000 

-23.14 

-0.37000 

-23.14 

^fcB  ~  fkB 

0.39948 

-32.17 

-0.39948 

-32.17 

FALSE 

0.41199 

-23.80 

-0.41199 

-23.80 

^fcA  ~  fkA 

0.45999 

-12.02 

-0.45999 

-12.02 

^fcA 

0.47000 

-13.75 

-0.47000 

-13.75 

VcB 

0.50000 

-21.61 

-0.50000 

-21.61 

^fcc  ~  fkc 

0.50800 

-21.79 

-0.50800 

-21.79 

^fcA  +  fkA 

0.57001 

-23.99 

-0.57001 

-23.99 

^fcB  d"  fkB 

0.59000 

-12.70 

-0.59000 

-12.70 

Voc 

0.67999 

-23.34 

-0.67999 

-23.34 

^fcc  d"  fkc 
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Table  15:  A+NS  (A=32768,  Ap=64),  15  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.04800 

-17.87 

-0.04800 

-17.87 

f^A 

0.09601 

-24.95 

-0.09601 

-24.95 

0.36401 

-27.89 

-0.36401 

-27.89 

'A 

(M 

1 

0.41199 

-22.47 

-0.41199 

-22.47 

2/cx  fkA 

0.45999 

-11.54 

-0.45999 

-11.54 

VcA 

0.50800 

-21.32 

-0.50800 

-21.32 

‘^fcA  d"  fkA 

0.55600 

-25.68 

-0.55600 

-25.68 

VcA  +  VkA 

Table  16:  B+N5  (A=32768,  Np=64),  9  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.10001 

-19.07 

-0.10001 

-19.07 

fkg 

0.37000 

-21.53 

-0.37000 

-21.53 

^fcB  ~  fkB 

0.47000 

-12.15 

-0.47000 

-12.15 

VoB 

0.57001 

-21.65 

-0.57001 

-21.65 

^fcB  d"  fkB 

A-7 


APPENDIX 


Table  17:  C+N5  (iV=32768,  Np=6i),  9  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.09000 

-18.75 

-0.09000 

-18.75 

fkc 

0.50000 

-21.52 

-0.50000 

-21.52 

^fcc  ~  fkc 

0.59000 

-12.17 

-0.59000 

-12.17 

Vcc 

0.67999 

-21.54 

-0.67999 

-21.54 

^fcc  +  fkc 

Table  18:  A+B+C+N5  (A=32768,  Ap=64),  27  Features  Detected. 


a 

Value  (dB) 

a 

Value  (dB) 

Feature 

0.00000 

0.00 

0 

0.04800 

-21.62 

-0.04800 

-21.62 

fkA 

0.09000 

-21.70 

-0.09000 

-21.70 

fkc 

0.09601 

-26.59 

-0.09601 

-26.59 

VkA 

0.10001 

-23.22 

-0.10001 

-23.22 

fkB 

0.37000 

-24.94 

-0.37000 

-24.94 

^fcB  ~  fkB 

0.41199 

-26.79 

-0.41199 

-26.79 

^fcA  ~  fkA 

0.45999 

-14.57 

-0.45999 

-14.57 

^fcA 

0.47000 

-16.26 

-0.47000 

-16.26 

VcB 

0.50000 

-24.36 

-0.50000 

-24.36 

^fcc  -  fkc 

0.50800 

-24.55 

-0.50800 

-24.55 

^fcA  d"  fkA 

0.57001 

-26.13 

-0.57001 

-26.13 

“^fcB  d"  fkB 

0.59000 

-16.03 

-0.59000 

-16.03 

Vcc 

0.67999 

-25.18 

-0.67999 

-25.18 

^fcc  d-  fkc 
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